First-Principles Perspective on Gas Adsorption by [Fe4S4]-Based Metal–Organic Frameworks

[Fe4S4] or [4S–4Fe] clusters are responsible for storing and transferring electrons in key cellular processes and interact with their microenvironment to modulate their oxidation and magnetic states. Therefore, these clusters are ideal for the metal node of chemically and electromagnetically tunable metal–organic frameworks (MOFs). To examine the adsorption-based applications of [Fe4S4]-based MOFs, we used density functional theory calculations and studied the adsorption of CO2, CH4, H2O, H2, N2, NO2, O2, and SO2 onto [Fe4S4]0, [Fe4S4]2+, and two 1D MOF models with the carboxylate and 1,4-benzenedithiolate organic linkers. Our reaction kinetics and thermodynamics results indicated that MOF formation promotes the oxidative and hydrolytic stability of the [Fe4S4] clusters but decreases their adsorption efficiency. Our study suggests the potential industrial applications of these [Fe4S4]-based MOFs because of their limited capacity to adsorb CO2, CH4, H2O, H2, N2, O2, and SO2 and high selectivity for NO2 adsorption.


Section S1. Computational details
The Gaussian 16 A.03 quantum chemical package 1 was used to study the thermodynamics and kinetics aspects of the studied [Fe4S4] clusters and metal-organic framework (MOF) models.
Since the Fe atoms in [Fe4S4] clusters can be found at different oxidative states and that the cluster might gain an antiferromagnetic state, the quantum chemical treatment of the systems was challenging. 2,3 To improve the accuracy of our results, the computational level was first calibrated using density functional theory (DFT)-based methods and basis sets to study the geometry of several experimentally resolved structures. It should be noted that, from the calibration step to the thermodynamics analysis and the reaction kinetics studies, all potential spin states were evaluated and the one with the minimum zero-point energy (ZPE) electronic energy (and Gibbs free energy) was chosen as the ground state. For computational level calibration, the ground spin state was evaluated at the PBE/def2-TZVP level and then the structures were optimized at various computational levels. For further thermodynamics and kinetics studies, the optimal computational level was used to identify the ground spin state.
Also, all calculations included the D3 Grimme's dispersion correction 4 because of the importance of dispersion effects for correct description of the reaction paths and reaction intermediates during the kinetics studies 3 and the potentially weak adsorbate-adsorbent interactions in the thermodynamics evaluations. The only exception was related to the ωB97XD-based calculations, which inherently includes dispersion impacts. 5 The first compound used in method calibration was [Fe4(NO)4(µ-S)4] -1 . As Table S1 reports, among the selected levels, M06-2X/def2-TZVP reproduces the experimental structure with the highest error (mean absolute error: MAE > 5%). However, PBE/def2-TZVP and PBE/6-311++G** give a significantly more accurate structure with only 0.7% mean absolute error (MAE). Therefore, only these two computational levels were applied to the next test molecules. Tables S2-S4 indicates that the two PBE/def2-TZVP and PBE/6-311++G** levels are equally efficient in reproducing the experimental geometries. However, PBE/6-311++G** (average MAE: 1.24%) slightly outperforms PBE/def2-TZVP (average MAE: 1.27%), also featuring a shorter CPU time (see Tables S1 and S2). Therefore, PBE/6-311++G** was selected as the optimal level. This computational level combines the choice of the PBE functional by Amitouche et al. 6 for studying gas adsorption on small Fe-S clusters and improves the choice of B3LYP/6-31++G** by Niu and Ichiye, 7 who stated that the addition of sp-type diffuse functions to the 6-31G** basis set enhances the accuracy of the redox energies calculated for [Fe4S4] clusters. Notably, S3 the success of the PBE functional in predicting the FeS cluster properties can be attributed to the cancellation of errors. To concern correlation effect in PBE-based calculations, one can use the DFT + Ueff scheme, in which Ueff is an orbital-dependent correction to the PBE functional or the Generalized Gradient Approximation (GGA). This parameter is generally expressed as the difference between the Hubbard U parameter, which is the Coulomb energetic cost to place two electrons at the same site, and an approximation of the Hund's exchange parameter J (i.e., Ueff = U -J). When U ≈ J, the correlation effect can be ignored and the GGA approach would be adequate. In the case of FeS clusters, Ueff is quite small. For example, the Ueff of Fe3S4 is about 1 eV, 8 and the U and J parameters of pure metallic Fe almost cancel completely. 9 Therefore, the cancellation of the U and J parameters on the iron sites leads to the significant performance of PBE for our studied system.   As an example, 13tet CLP was found 59.5 kJ mol -1 more stable than its singlet antiferromagnetic state, at 0 K. Furthermore, the difference in their Gibbs free energy at 298.15 K and 1 atm was calculated to be 82.5 kJ mol -1 . Therefore, we skipped the antiferromagnetic state in our thermodynamics and kinetics studies, but still screened all potential ferromagnetic spin states when the presence of the adsorbate could alter the spin state of the system (such as in 3 O2 and 2 NO2 adsorption).

Comparison of the optimized geometries with the experimental structures in
After calibrating all basic quantum chemical parameters, we started the thermodynamics and In addition to the thermodynamics analysis, we estimated the oxidative and hydrolytic stability of the clusters and MOF models by studying their reaction with 3 O2 and 1 H2O, respectively.
First, a number of different reaction transition states (TSs) were guessed and optimized. If the optimized structure associated with one imaginary frequency, intrinsic reaction coordinate (IRC) analysis 22,23 was performed to connect that TS to the related product(s) and reactant(s).
Often, the proposed TSs did not lead to any product/reactant and were neglected. The identified reactants/products were optimized and accepted as reactants, products, or reaction S7 intermediates if they lacked any imaginary frequencies. Then, the confirmed reaction paths were combined to construct the potential energy surfaces (PESs).  -105.2, -98.9 -177.5, -174.5 -51.7 to -47.2 -50.5 to -25.8 a The min-max value ranges reported indicate several unique adsorption modes. Similarly, single or two discrete values indicate that the geometry or energy of two or several starting adsorption configurations have converged to the same geometry/energy value.  Figure S1. The lowest energy (most favorable) adsorbate/adsorbent configurations.